#include "findMisMatch.h"



__global__ void baseMisMatch(char *d_read, char *d_ref, unsigned int *d_readPos, unsigned int *d_refPos, unsigned char* d_numMisMatch, unsigned char** d_misPos, unsigned char** d_misBase, const int readLength)
{
        __shared__ unsigned int misMatch[FQUERY_PER_BLOCK];
	unsigned int x = threadIdx.x;
	unsigned int i = blockIdx.x*FQUERY_PER_BLOCK  + x;
	unsigned int offset = threadIdx.y;
	unsigned int refPos = d_refPos[i];
	unsigned int readPos = d_readPos[i];
	misMatch[x] = 0;
	while(offset < readLength)
	{      
		if(d_ref[refPos+offset] != d_read[readPos+offset])
	        {      
		   unsigned int old =  atomicInc(&misMatch[x], FWARP_SIZE)%2;
		   d_misPos[old][i] = offset;
		   d_misBase[old][i] = d_read[readPos+offset];
		     
	        }
		offset += FWARP_SIZE;
        }
	d_numMisMatch[i] = misMatch[x];

}
__global__ void multiplyLength(unsigned int *d_readPos, const int readLength)
{
	d_readPos[blockIdx.x*FQUERY_PER_BLOCK + threadIdx.x] *= readLength;
}
void handleError(cudaError_t err)
{
	if(err != cudaSuccess)
	{
		printf("error occur: %s\n", cudaGetErrorString(err));
	}
}
void findMisMatch(char* read, char* ref, const int readLength, unsigned int* queryId, unsigned int *pos, unsigned int numElement, struct Info* info)
{
	unsigned char **d_misPos;
	unsigned char **d_misBase;
	unsigned char *tempBase[2];
	unsigned char *tempPos[2];
	unsigned char *numMisMatch = info->numMisMatch;
	unsigned char *misBase1 = info->misBase[0];
	unsigned char *misBase2 = info->misBase[1];
	unsigned char *misPos1 = info->misPos[0];
	unsigned char *misPos2 = info->misPos[1];
	char *d_read;
	char *d_ref;
	unsigned int *d_readPos;
	unsigned int *d_refPos;
	unsigned int readLen = info->readStringLength;
	unsigned int refLen = info->refStringLength;
	unsigned int dataBlock = FBLOCK_NUM*FQUERY_PER_BLOCK;
	unsigned char* d_numMisMatch;
	unsigned char* d_misBase1;
	unsigned char* d_misBase2;
	unsigned char* d_misPos1;
	unsigned char* d_misPos2;
	size_t allocateSize;
	size_t uintSize = sizeof(unsigned int);
	size_t biCharPoints = sizeof(unsigned char*)*2;
	float calTime = 0;
	float sumTime = 0;
	float memTime = 0;

	unsigned int alloc[2];
	unsigned int pass = numElement/dataBlock;
	alloc[1] = numElement - pass*dataBlock;
	alloc[0] = dataBlock;
	dim3 blockSize(FQUERY_PER_BLOCK, FWARP_SIZE);
	
	timeval memStart, memEnd;
	timeval calStart, calEnd;
	timeval sumStart, sumEnd;
	gettimeofday(&memStart, NULL);
	gettimeofday(&sumStart, NULL);

	printf("findMisMatch: processing...\n");

	cudaFuncSetCacheConfig(baseMisMatch, cudaFuncCachePreferL1);
	cudaFuncSetCacheConfig(multiplyLength, cudaFuncCachePreferL1);
	
	cudaMalloc(&d_read, readLen+refLen+1);
	d_ref = d_read+readLen;

	cudaMalloc(&d_readPos, dataBlock*uintSize*2);
	d_refPos = d_readPos + dataBlock;
	
	cudaMalloc(&(d_numMisMatch), dataBlock*5);
	d_misBase1 = d_numMisMatch + dataBlock;
	d_misBase2 = d_misBase1 + dataBlock;
	d_misPos1 = d_misBase2 + dataBlock;
	d_misPos2 = d_misPos1 + dataBlock;
	
	cudaMalloc(&d_misPos, biCharPoints);
	cudaMalloc(&d_misBase, biCharPoints);

	tempBase[0] = d_misBase1;
	tempBase[1] = d_misBase2;
	tempPos[0] = d_misPos1;
	tempPos[1] = d_misPos2;

	cudaMemcpy(d_misPos, &tempPos, biCharPoints, cudaMemcpyHostToDevice);
	cudaMemcpy(d_misBase, &tempBase, biCharPoints, cudaMemcpyHostToDevice);
	cudaMemcpy(d_read, read, readLen, cudaMemcpyHostToDevice);
	cudaMemcpy(d_ref, ref, refLen, cudaMemcpyHostToDevice);

	cudaThreadSynchronize();
	
	gettimeofday(&memEnd, NULL);
	gettimeofday(&calStart, NULL);
	

	for(unsigned int i = 0; i <= pass; i++)
	{
		allocateSize = alloc[i/pass];
		cudaMemcpy(d_readPos, queryId + i*dataBlock, allocateSize*uintSize, cudaMemcpyHostToDevice);
		cudaMemcpy(d_refPos, pos + i*dataBlock, allocateSize*uintSize, cudaMemcpyHostToDevice);
		
		multiplyLength<<<FBLOCK_NUM, FQUERY_PER_BLOCK >>>(d_readPos, readLength);
		cudaThreadSynchronize();

		baseMisMatch<<<FBLOCK_NUM, blockSize>>>(d_read, d_ref, d_readPos, d_refPos, d_numMisMatch, d_misPos, d_misBase, readLength);	
		
		cudaMemcpy(numMisMatch+i*dataBlock, d_numMisMatch, allocateSize, cudaMemcpyDeviceToHost);
		cudaMemcpy(misPos1+i*dataBlock, d_misPos1, allocateSize, cudaMemcpyDeviceToHost);
		cudaMemcpy(misPos2+i*dataBlock, d_misPos2, allocateSize, cudaMemcpyDeviceToHost);
		cudaMemcpy(misBase1+i*dataBlock, d_misBase1, allocateSize, cudaMemcpyDeviceToHost);
		cudaMemcpy(misBase2+i*dataBlock, d_misBase2, allocateSize, cudaMemcpyDeviceToHost);

	}
	
	cudaThreadSynchronize();
	gettimeofday(&calEnd, NULL);

	cudaFree(d_read);
	cudaFree(d_ref);
	cudaFree(d_readPos);
	cudaFree(d_refPos);
	cudaFree(d_numMisMatch);
	cudaFree(d_misBase1);
	cudaFree(d_misBase2);
	cudaFree(d_misPos1);
	cudaFree(d_misPos2);

	cudaThreadSynchronize();
	
	gettimeofday(&sumEnd, NULL);

	calTime = (float)(FCLK*(calEnd.tv_sec - calStart.tv_sec)+(calEnd.tv_usec - calStart.tv_usec))/FCLK;
	sumTime = (float)(FCLK*(sumEnd.tv_sec - sumStart.tv_sec)+(sumEnd.tv_usec - sumStart.tv_usec))/FCLK;
	memTime = (float)(FCLK*(memEnd.tv_sec - memStart.tv_sec)+(memEnd.tv_usec - memStart.tv_usec))/FCLK;
	printf("findMisMatch: finished\n");
	printf("findMisMatch: Total time: %.2f s\n", sumTime);
	printf("findMisMatch: Time used copying data from Host to Device: %.2f s\n", memTime);
	printf("findMisMatch: Calculation time: %.2f s\n", calTime);
}
